classdef Chebyshev2 < FilterFather
    properties
    end
    methods
        
        function this = Chebyshev2(guiData)
            this.calculateFilter(guiData);
            minEdge = min([guiData.Fp guiData.Fa])*2*pi;
            maxEdge = max([guiData.Fp guiData.Fa])*2*pi;
            this.calcData([minEdge maxEdge]);
        end

        function calculateFilter(obj,guiData)
            Aa = guiData.Aa;
            Ap = guiData.Ap;
            obj.Approximation = 'Chebyshev II';
            obj.FilterType = guiData.filterType;
            Wp = guiData.Fp*2*pi;
            Wa = guiData.Fa*2*pi;
            [Wp Wa] = obj.adjustBand(Wp,Wa);
            obj.DesnormValue = guiData.DenormVal;
            obj.DesnormType = guiData.DenormType;
            fixOrder = guiData.fixOrder;
            MaxOrder = guiData.order;
%             MaxQ = guiData.maxQ;
            obj.normalization(Wa,Wp);
%             obj.Wdesnor = Wp;
            
            
            epsilon = 1/sqrt(10^(.1*abs(Aa))-1);

            % Order calculus
            N = ceil(acosh( 1/(sqrt(10^(.1*abs(Ap))-1)*epsilon) ) /acosh(1/obj.WaNormalized));
            
            if MaxOrder == 0
                obj.Order = N;
            else 
                if(N > MaxOrder || fixOrder == true)
                    N = MaxOrder;
                end
                obj.Order = N;
            end
%             if obj.Order > MaxOrder
%                 obj.Order = MaxOrder;
%             end
            N = obj.Order;
            
            a = sinh(1/N*asinh(1/epsilon));
            b = cosh(1/N*asinh(1/epsilon));
            %x = a*sin(pi*(1:2:N-1)/(2*N));
            %y = b*cos(pi*(1:2:N-1)/(2*N));
            %p = (x - 1j*y)./(x.^2+y.^2);
            x = -1*a*sin(pi*(1:2:N-1)/(2*N));
            y = b*cos(pi*(1:2:N-1)/(2*N));
            p = 1./(x + 1i.*y);
            p = [p, conj(p)]; % To sort numerical errors
            p = p(:);
            if rem(N,2)==1   % N is odd
                p = [p; -1/a];    % The negative real pole
            end
%             K = epsilon*cosh(N*acosh(obj.WaNormalized));

%             [x y] = pol2cart(atan(-cot(pi*(2*(1:1:N)-1)/(2*N))*coth(asinh(K)/N))+pi,obj.WaNormalized./sqrt((cos(pi*(2*(1:1:N)-1)/(2*N))).^2+(sinh(asinh(1:1:N)/N)).^2));
%             p = x+1i*y; % To sort numerical errors
%             
            z=1j./(cos(pi*((1:2:N-1))/(2*N)));
            z = [z conj(z)];
            % si N es impar hay un cero en infinito 
            

            obj.Poles = p;
            obj.Zeros = z;
            obj.K = abs(real(prod(obj.Poles)/prod(imag(obj.Zeros))));
            obj.frecDesnormInvCheb(Ap);
            obj.K = abs(real(prod(obj.Poles)/prod(obj.Zeros)));
            obj.desnormalization();
            obj.TransferFunction = obj.TransferFunction*10^(guiData.gain/20);
        end
        
        function frecDesnormInvCheb(obj,Ap)
            num = poly(obj.Zeros);
            den = poly(obj.Poles);
            H = tf(num,den);
            [mag] = bode(H,1e-6);
            H = H/mag;
            [mag,~,w] = bode(H,linspace(0.1*obj.WaNormalized,2,10000));
            mag = -20*log10(mag(:,:));
            k = find(mag < Ap,1,'last');
            Wadj = (obj.WaNormalized-w(k))/w(k);
            obj.DesnormValue = 100 - obj.DesnormValue;
%             
            switch obj.DesnormType
                case 'Fp'                   
                    obj.Zeros = obj.Zeros*(1 + Wadj);
                    obj.Poles = obj.Poles*(1 + Wadj); 
                case 'Fother'
                    Wadj = Wadj*obj.DesnormValue/100;
                    obj.Zeros = obj.Zeros*(1 + Wadj);
                    obj.Poles = obj.Poles*(1 + Wadj);    
            end
        end
    end
end